Flow analysis from multiparticle azimuthal correlations 
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We present a new method for analyzing directed and elliptic flow in heavy ion collisions. Unlike 
standard methods, it separates the contribution of flow to azimuthal correlations from contribu- 
tions due to other effects. The separation relies on a cumulant expansion of multiparticle azimuthal 
correlations, and includes corrections for detector inefficiencies. This new method allows the mea- 
surement of the flow of identified particles in narrow phase-space regions, and can be used in every 
regime, from intermediate to ultrarelativistic energies. 
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_ ■ I. INTRODUCTION 

Qh' 

In noncentral heavy ion collisions, it is possible to measure azimuthal distributions of outgoing particles with 
ly— ^ . respect to the reaction plane. This is the so-called "flow analysis" , which is being actively studied over a wide range 
' of colliding energies, from below 25 MeV per nucleon in the center-of-mass system Q to above 60 GeV (|]. The main 
motivation for such studies is that anisotropics in the azimuthal distributions are likely to contain much information 
CN| ' on the physics in the hot, dense central region of the collision (sec for reviews). In particular, they may provide 
a signature of the formation of a quark-gluon plasma at ultrarelativistic energies Azimuthal distributions may 

also be of interest in connection with more exotic phenomena, such as the formation of disoriented chiral condensates 
IIP], or the study of parity and/or time-reversal violation |lC[ ], Finally, the combination of flow analysis and two- 
particle interferometry yields a three-dimensional picture of the emitting source |l3[| , as was shown recently at the 
' Brookhaven Alternating Gradient Synchrotron (AGS) H]. 

In this paper, we propose a new method to measure azimuthal distributions. As usual, they will be characterized 
' by their Fourier coefficients Il5|] 
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v n = ( e *"W-**)) = (cosn(0- (1) 



where </> is the azimuthal angle of an outgoing particle measured in the laboratory coordinate system, is the 



azimuth of the impact parameter (or reaction plane) and angular brackets denote a statistical average, over many 
events. The first two coefficients v\ and Vi are usually referred to as directed flow and elliptic flow, respectively. 
The u„'s can be measured for various particle species, as a function of transverse momentum and/or rapidity: we 
refer to these detailed measurements as to "differential" flow, following fusl - In this paper, we also discuss global 
measurements of v n , averaged over a large phase space region, typically corresponding to the acceptance of a detector. 
d ' We call this "integrated" flow. 

Since the reaction plane cannot be measured directly, the only way to obtain the coefficients v n experimentally is 
to deduce them from the azimuthal correlations between the outgoing particles: the correlation between every particle 
and the reaction plane induces correlations among the particles (which we call hereafter the "flow correlations"), from 
which v n can be reconstructed. The method we propose here is based on a systematic analysis of multiparticle 
azimuthal correlations. 



The most widely used method for the flow analysis is that initially proposed by Danielewicz and Odyniec |17| (see 
also for further developments), which relies on azimuthal correlations between two "subevents" . It has recently 

been applied at intermediate energies in Darmstadt Jl],|2l| , and at higher energies in Dubna p2]-|2^] , at the Brookhaven 
AGS P5H27|, at the CERN Super Proton Synchrotron (SPS) [§8|,|9| and, finally, at the Brookhaven Relativistic Heavy 
Ion Collider (RHIC) ||. An alternative, simpler method extracts flow from two-particle correlations and is still 
in use, both at intermediate |3l]] and at ultrarelativistic energies |32|,|33|]. Both methods are more or less equivalent: 
correlating two subevents amounts to summing two-particle correlations. In these analyses, one usually assumes that 
the only sources of azimuthal correlations are flow and, when necessary, transverse momentum conservation |34j]. 

However, we have shown in two papers p5| , |36"| that this assumption is no longer valid at SPS energies, where 
"direct" , nonflow two-particle correlations become of the same magnitude as the correlations due to flow. Even when 
nonflow correlations are a priori smaller than flow correlations, they must be taken into account in order to obtain 
accurate and reliable results. Some sources of nonflow correlations are well known. One can attempt to avoid them 
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experimentally by appropriate cuts in phase s pace M , or one can take them into account in the analysis , as was 
done for transverse momentum conservation |3^j36t] , for correlations from tto — > 77 decays f37jj38| , and for quantum 
correlations [ p5[ . But there is no systematic way to separate the effects of flow from other effects at the level of 
two-particle correlations. 

There have been several attempts in the past to go beyond two-particle correlations: analyses of the global event 
shape (3^] allowed the first observations of collective flow at intermediate 40 and ultrarelativistic energies |4l|]42"f , 
which were not biased by nonflow correlations. Multiparticle azimuthal correlations were also used in p3j. These 
methods are now considered obsolete because they apply only to the "integrated" flow, as defined above, whereas 
most of the recent analyses concentrate on the differential flow of identified particles, in particular kaons J26p4] ], 77 
mesons f21|| , A hyperons 27,4f|, and antiprotons 

In recent papers , we have shown for the first time that nonflow correlations can be removed systematically not 

only for integrated flow, but also in analyses of differential flow. However, a limitation of this method when measuring 
v n is the interference with higher harmonics (i>2n, v 3n, etc.). This interference may hinder the measurement of directed 
flow when elliptic flow is larger, which is likely to be the case at RHIC energies p9[ . 

Here we present an improvement of this method, which is free from this limitation, and in many respects simpler. In 
particular, it no longer involves the event flow vector on which most analyses are based As m our previous 

method, we perform a cumulant expansion of multiparticle azimuthal correlations, which eliminates order by order 
nonflow correlations, and can be used even if the detector does not have full azimuthal coverage. 

In Sec. ||, we show how to construct the cumulants of multiparticle azimuthal correlations by means of a generating 
function. These cumulants allow us to reconstruct the integrated flow from the measured correlations. The method 
is extended to differential flow in Sec. III. The relation with other methods is discussed in Sec. [V. Results of 



Monte-Carlo simulations are presented in Sec. V|. The most technical points are left to appendices: the construction 
of cumulants is explained in detail in Appendix [A|; interpolation formulas used to obtain the cumulants from the 
generating function are given in Appendix |Bj acceptance corrections, which extend the validity of the method to 
detectors with partial azimuthal coverage, are derived in Appendix finally, statistical errors on the flow values 
deduced from the cumulants are evaluated in Appendix 

The essential improvement on our previous method is the use of a new generating function, defined in Sec. II B, 
which corrects the limitations encountered in These improvements are discussed in detail in Sec. IV and in 

Appendix |A]; they are seen clearly in the simulations presented in Sec, 
acceptance corrections (Appendix and statistical errors (Appendix 

apply to our previous method. Apart from these differences, most of the material discussed in Sees, [n] and III can be 
found in |47]|, although the present derivation is more transparent. 



In addition, the detailed discussions of 
are completely new, although they also 



II. INTEGRATED FLOW 



In Sec. II A, we illustrate with a few examples the principle of the cumulant expansion of multiparticle azimuthal 
correlations, and show how it can be u sed to perform flow measurements with a better sensitivity than the previous 
metho ds. T hen we explain, in Sec. II B , how to perform this expansion in practice, by means of a generating function. 
In Sec. IIC, we derive the relations between the cumulants and the flow v n , integrated over some phase space region. 
Using cumulants to various orders, one thus obtains different estimates for v n . The uncertainties associated with each 
estimate due to nonflow co rrela tions and limited statistics, and the resulting optimal choice, are examined in Sec. II D. 
Finally, we discuss in Sec. HE the generalization of the previous subsections to different, optimal particle weights. 



A. Cumulants of multiparticle azimuthal correlations 



We denote by tfij, with j — I,...,M, the azimuthal angles of the particles seen in an event with multiplic- 
ity M, measured with respect to a fixed direction in the detector (this was denoted by <f>j in |47j]). In this 
paper, we shall be concerned with multiparticle azimuthal correlations, which we write generally in the form 
(exp[m(0i + ■ ■ ■ + 4>k — 4>k+i — ■ ■ ■ — where n is the Fourier harmonic under study (n = 1 for directed flow, 

n = 2 for elliptic flow), and the brackets indicate an average which is performed in two steps: first, one averages over 
all possible combinations of k + 1 particles detected in the same event; then, one averages over all events. 

Correlations between k + 1 particles can be generally decomposed into a sum of terms involving correlations between 
a smaller number of particles. Consider for instance the measured two-particle azimuthal correlation <V"W 1- ^ 2 )). It 
can be written as 
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^ e in(0i-02)\ — {V"^ 1 ) (e~" l< ^ 2 ) + ^e™^ 1- ^ 2 -'^} (2) 

where le'™^ 1 "^ 2 ') is by definition the second order cumulant. To understand the physical meaning of this quantity, 
we first consider a detector whose acceptance is isotropic, i.e., which does not depend on 4>. Such a detector will be 
called a "perfect" detector. Then, the average (e m ^A vanishes by symmetry [since <pj is measured in the laboratory, 
not with respect to the reaction plane, ^e m ^ J ) does not correspond to the flow v n defined in Eq. (jjj)]: the first term on 
the right-hand side (r.-h. s.) of Eq. (EJ) vanishes and the cumulant reduces to the measured two-particle correlation. 

The relevance of the cumulant appears when considering the more realistic case of a detector with uneven acceptance. 
Then, the first term on the r.-h. s. of Eq. (Eft can be nonvanishing. But the cumulant vanishes if <fii and <f>2 
are uncorrelated. Thus the cumulant ^e 1 ™^ 1- ^ 2 )} isolates the physical correlation, and disentangles it from trivial 
detector effects. 

There are several physical contributions to the correlation ^e m ^ 1_< ^ 2 ^)), which separate into flow and nonflow 
(or direct) correlations. When the source is isotropic (no flow), only direct correlations remain. They scale with the 
multiplicity M like 1 /M J35],[36| , as can be easily understood when considering correlations between the decay products 
of a resonance: when a p meson decays into two pions, momentum conservation induces an angular correlation of 
order unity between the decay pions; besides, the probability that two arbitrary pions seen in the detector result from 
the same p decay scales with the total number of pions like 1/M. All in all, the correlation between two arbitrary 
pions is of order 1/M . If the source is not isotropic, flow, which is by definition a correlation between emitted particles 
and the reaction plane, generates azimuthal correlations betw een a ny two outgoing particles, and gives a contribution 



d„ to the second-order cumulant, as will be explained in Sec. [IC. One can measure the flow using the second-order 
cumulant if this contribution dominates over the nonflow contribution, i.e., if v n » \j\[M |35||36|. This is the domain 
of validity of standard flow analyses, which are based on two-particle correlations. 

Our main point is that through the construction of higher order cumulants, one can separate flow and nonflow 
correlations. To illustrate how this works, we consider for simplicity a perfect detector. Then, we decompose the 
measured four-particle correlation as follows: 



e -.. (01 +02-03-04) \ _ /g"i(0l-03)\ / e «n(02-04)\ _|_ /gj«(0l-04)\ /gm(0 2 -0 3 )\ 4_ // p t"(01+02 



(3) 



If the particles are correlated pairwise, there are two possible combinations leading to a nonvanishing value for the 
left-hand side (l.-h. s.): the pairs can be either (1,3) and (2,4), or (1,4) and (2,3). This yields the first two terms 
in the right-hand side. The remaining term |e m ^ 1+ ^ 2_ ^ 3_c!i4 ^)), which is by definition the fourth-order cumulant, 
is thus insensitive to two-particle nonflow correlations. However, it may still be influenced by higher order nonflow 
correlations: if, for instance, a resonance decays into four particles, the resulting correlations between the reaction 
products do not factorize as in Eq. (||). We call such correlations "direct" four-particle correlations. Fortunately, 
their contribution to the fourth-order cumulant is very small: it scales with the multiplicity like 1/M 3 |^7|, while the 
measured corr elation ( e i ™Wi+-fe-0 3 -04)^ [ s generally much larger, of order 1/M 2 [the two -particle correlation terms in 
the r.-h. s. of Eq. (||) are of orde r 1/M , as explained above]. On the other hand, flow yields a contribution — to the 



cumulant, as we shall see in Sec. II C . Therefore, the cumulant is dominated by the flow as soon as v n ^> 1 /A/ 3 / 4 . This 
is a major improvement on two-particle correlations, which are limited by the much stronger constraint v n 3> 1/VM. 
Equation (^) can be rewritten 

^gi«(01+02-03-04)^ _ 2 ^g™(01-03)^ 2 _|_ ^g'in(01 +02-03-04)^ ^ 

where we have used the symmetry between <pi and 4>2 (resp. 03 and ^4). However, Eqs. (^|) and (||) only hold for a 
perfect detector, therefore they are of little practical use. It is in fact possible to build an expression for the fourth- 
order cumulant which eliminates both detector effects and nonflow correlations, but this expression is very lengthy. 



This is the reason why we introduce a generating function of cumulants in Sec. II B. It will enable us to construct 
easily cumulants of arbitrary orders for arbitrary detectors. 

More generally, the cumulant ((exp[m(<^i + • • • + <pk — <fik+i — ■ ■ • — (f>k+i)]}) , which involves k + 1 particles, is of order 
M 1 ~ l when there is no flow. It eliminates all nonflow correlations up to order k + I — 1. Only direct correlations 
between k + 1 particles remain. Cumulants with k =/= I vanish for a perfect detector and are physically irrelevant. The 
interesting cumulants are the "diagonal" ones, with k — I, as in E qs. ( %) and (|J). The contribution of flow to these 

When this contribution dominates over the 



cumulants, proportional to z; 2fc , will be evaluated precisely in Sec. |IC. 
nonflow contribution, the measured cumulant yields an estimate of the value of v n , which we denote by v n {2k}, where 
k > is in principle arbitrary. 
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B. Generating function 

Cumulants can be expressed elegantly, and without assuming a perfect detector as in Eq. (Q), using the formalism 
of generating functions. For each event, we define the real- valued function G n (z), which depends on the complex 
variable z — x + iy: 



M 



Gn{z) = n 1 



W 

J = l 

= TTM-| 2xcos ( n ^i) +2j/sin(n^) \ ^ 
j= i V M J 

where z* = x — iy denotes the complex conjugate. This generating function can then be averaged over events with the 
same multiplicity M. We denote this statistical average by (G n (z)). Its expansion in power series generates measured 
azimuthal correlations to all orders: 

(Gn(z)) =1 + ^:^ 




\j<k 

1 + z(e~ m ^) + z* (e m ^) + ^± ( — / e -*»(*i+fc) \ + 1_ / e *»(*+*0\ + zz * / e in(^-«\ j + ... (6) 



where, the averages (e™^ 1 ), (e" 1 ^ 1 ^ 2 '), etc. are defined as in Sec. II A. More generally, expanding (G n (z)) to order 
z* k z l yields, up to a numerical coefficient, the (£; + Z)-particle correlation (exp[m(0i + ••• + <j>k — <f>k+i — • • • — <pk+i)]). 
The generating function (G n (z)) thus contains all the information on measured multiparticle azimuthal correlations. 

If the detector is perfect, the statistical average (G n (z)) does not depend on the phase of z; it only depends on 
\z\ = y 1 x 2 + y 2 . To see this, one may note that changing z into ze™ in the generating function (|^) amounts to 
shifting all angles by the same quantity c/)j — * <j)j — a. Now, the probability for an event to occur is unchanged under 
a global rotation; therefore (G n (z)) is unchanged under this rotation, hence the result. In this case, the only terms 
which remain in the power-series expansion are the isotropic terms, proportional to z k z* k , which involve only relative 
angles. 

The generating function provides us with a way to obtain a compact expression for cumulants of arbitrary orders. 
We define the generating function of the cumulants C n (z) by 

C n (z) = M ((G„(z)) 1/M -l) . (7) 

The expansion of this function in power series of z and z* defines the cumulants: 

C n {z) = ((e in W 1+ -+**-**+ 1 --* k +'))) . (8) 

One easily checks that if the particles are uncorrelated, all the cumulants vanish beyond order one, i.e., for k + 1 > 2. 
Indeed, if all the m Eq. (||) are independent from each other, the mean value of the product is the product of the 
mean values, so that 

The generating function of cumulants, Eq. (Q), then reduces to 

C n (z) = z*(e m ^) +z(e- m ^). (10) 
Comparing with Eq. (0), cumulants of order 2 and higher vanish when particles are uncorrelated, as expected. 
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The cumulant ^e 1 "^ 1 obtained when expanding Eqs. (^) and (||) to order zz* coincides with the second-order 

cumulant defined in Eq. (||) in the limit of large M (see Appendix [A|). Expanding C n (z) to order z 2 z* 2 , one obtains an 
expression for the cumulant ^ e ™(.<l>i+<fa-fo-<t>*)^ which reduces to Eq. (Q) for a perfect detector. But the expression 
derived from Eqs. (]?]) a nd (p[ ) is still valid with an imperfect detector, while Eq. (Q) is not. 

As mentioned in Sec. pAT cumulants with k ^ I vanish for a perfect detector, since the generating function C n (z) 
in Eq. (||) depends only on \z\. The interesting cumulants are the diagonal terms with k = I, which are related to the 
flow. We denote them by c n {2k}: 

c n {2k} ee ^e m(01+ -+^~^+ 1 02fc) Jj> . (11) 

In practice, expanding the generating function C n {z) analytically is rather tedious beyond order 2. The simplest way 
to extract c„{2fc} is to tabulate the generating function (Q), and then compute nu meri cally the coefficients of its 



power-series expansion, using interpolation formulas which can be found in Appendix B 1 . 

Finally, we have assumed here that the multiplicity M is exactly the same for all events involved in the analysis. In 
practice, one performs the flow analysis for a class of events belonging to the same contrality interval, and M fluctuates 
from one event to the other. That explains our introducing the factor 1/M in the definition of the generating function 
(|J), as explained in more detail in Appendix |a|. The average over events (G n (z)) then involves an average over M, 
and M must be replaced by its average value (M) in the definition of the cumulants, Eq. (0). This, however, leads to 
errors, especially when the acceptance is bad (see Appendix This point will also be illustrated by the simulations 
presented in Sec. 0. In order to avoid these effects, one may use only a randomly selected subset of the detected 
particles, with a fixed multiplicity M, to construct the generating function. 



C. Contribution of flow to the cumulants 



Let us evaluate the contribution of flow to the cumulants c n {2k}. For sim plici ty, we assume that the detector is 
perfect. The generalization to an uneven acceptance is performed in Appendix |C l| . Under this assumption, one easily 
calculates the generating function (G n (z)) and, from it, the values of the cumulants. 

Let us call §r the azimuthal angle of the reaction plane of a given event. The average over events can be performed 
in two steps: one first estimates the average over all events with the same then one averages over Qr. We denote 
by (x|$ij) the average of a quantity x for fixed With this notation, the definition of v n , Eq. (|l]), gives 

{e in ^\^ R ) =v n e in * R . (12) 

Neglecting for simplicity nonflow correlations between particles, we obtain from Eq. (^]) 

z v n e ~ m *« + z*v n e™* B N M 



(G n (z)\<f R ) = [l + — ) . (13) 



One must then average over 



'^O = I*" {G n (z)\$ R )^f. (14) 



Inserting Eq. ( |T3| ) in this expression, one obtains 



(G n (z)) - g (M _ 2fe)!(jfc!)2 © N 2 * 

^I (2v n \z\), (15) 



where, in the last equation, we have assumed that M is large, so that M\/{M — 2k)\ = M 2k and one may extend the 
sum over k to infinity. Iq denotes the modified Bessel function of the first kind. The result depends only on \z\, as 
expected from the discussion in Sec. II C . 

The generating function of the cumulants (j?]) now reads 



C 



n (z) ~ M (/ (2igz|) 1 / M - l) ~ ln/ (2f„ \z\). (16) 
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This equation can be expanded in power series. Comparing with Eq. (g), the cumulants with k ^ I vanish, as expected 
for a perfect detector, while the diagonal cumulants c n {2k} defined by Eq. are related to v n . From the measured 
c„{2fc}, one thus obtains an estimate of v n , which is denoted by v n {2k}. The lowest order estimates are 

v n {2} 2 ee c„{2}, (17a) 
t>„{4} 4 ee -c„{4}, (17b) 
v n {6} e ee c„{6}/4. (17c) 

When the detector acceptance is far from isotropic, as is the case of the PHENIX detector at RHIC |33|], which covers 
approximately half the range in azimuth, these relations no longer hold. The issue of acceptance corrections, discussed 
in detail in Appendix |c| is more subtle than might be thought at first sight, for the following reason: when there is 
flow, the probability that a particle be detected depends on the orientation of the reaction plane &r if the detector 
only has partial azimuthal coverage. Hence, if a fixed number of particles are emitted, the number of particles seen in 
the detector depends of Reciprocally, for a fixed value of the multiplicity M seen in the detector, the probability 
distribution of is not uniform, which creates an important bias in the flow analysis. In the calculations of Appendix 



C 1 , we assume the centrality selection is done by an independent detector (for instance a zero-degree calorimeter) 
which has (at least approximately) full azimuthal coverage, so that the distribution of $ij is uniform for the sample 
of events used in the flow analysis. 

Under this assumption, one can derive general relations between the cumulants and the flow. It turns out that, 
in general, c n {2k} depends not only on v n , but also on other harmonics v p with p ^ n. In order to obtain the 
corresponding relations, we introduce the acceptance function A(<p), which is the probability that a particle with 
azimuthal angle <j) be detected. The Fourier coefficients of this acceptance function are 

/■27T J / 



e-^A{4>)^. (18) 

The relati ons between the cumulants c n {2k} and the estimates v p {2k} involve these coefficients. They are deriv ed i n 
Appendix CI, to leading order in v p . The results for directed flow and elliptic flow are given by Eqs. ( |5^ ) and JC7|), 



respectively. 



D. Errors 



We now examine the orders of magnitude of systematic errors, arising from unknown nonflow correlations, and 
statistical errors, due to the finite number of events available. More precisely, we estimate the difference Sv n \2k} 
between the true integrated flow v n and its values reconstructed from the cumulants, v n {2k}, defined in Eqs. (fFfl). 
We show which value of 2k minimizes the total uncertainty. 



As explained in Sec. [I A, nonflow 2fc-particle correlations give a contrib ution of order M 1 2k to the cumulant 
c n {2k}. This is to be compared with the contribution of flow derived in Sec. IIC, of order v 2 ^. We may thus write 



v n {2k} 2k - v 2k = 0(M ), (19) 

which is an estimate of the systematic error 8v n {2k} due to nonflow correlations. Obviously, flow can be measured 
only if v 2 ^ ^S> M 1 ~ 2k , For large orders k ^> 1, this condition becomes 

v n > 1/M, (20) 

which is a necessary condition for the flow to be measurable [|^J . We believe there is no way to extract a flow of order 
1/M or smaller. 

In this paper, we always assume that condition (^0|) is fulfilled. If this is the case, the systematic error on v n given 
by Eq. (|l9|), Sv n {2k} ~ {Mv n ) 1 ~ 2k , becomes smaller and smaller as k increases: thus one should construct cumulants 
of orders as high as possible. 

One must also take into account the statistical error, due to the finite number of events N cvts available. The order 
of magnitude of statistical errors can easily be understood. The cumulant c n {2k} involves correlations between 2k 
particles belonging to the same event. There are roughly M 2k ways (for large enough M) to choose 2k particles 
among the M particles detected, and one averages over all possible combinations. Since this is done for all iV ev t s 
events, there is a total of M 2k N cvts subsets of 2k particles involved in the evaluation of the cumulants. The resulting 
statistical error is therefore 



G 



v n {2k} 2k -vf = 0\ . (21) 

Unlike the systematic error, the statistical error genera lly i ncreases with increasing cumulant order 2k (it may in fact 



decrease in some cases, but only slightly, see Appendix D 2). Therefore, the order 2k which gives the best compromise 



is that for which statistical and systematic errors are both of the same magnitude. Equating the right-hand sides of 
Eqs. @ and (||), one obtains the optimal cumulant order Jl7|: 

^ H3 + 1F- (22) 

In most practical cases, the fourth-order cumulant (2k = 4), that is, removing two-particle nonflow correlations, is to 
be preferred. 



Statistical errors are discussed more thoroughly in Appendix D 2. There, we derive exact formulas for the standard 
deviations of the cumulants, and for their mutual correlations. Two regimes can be distinguished, depending on the 
value of the dimensionless parameter \ = v n \^M 7 which has been used previously as a measure of the reaction plane 
resolution [^o[ . If x <C 1, the standard errors agree with the simple estimate (pl|), and different estimates v n {2k} and 
v n {2l} with k I are uncorrelated. If \ ^S> 1, on the other hand, they are strongly correlated and the standard error 
becomes 

S(v n {2k}) stat = 1 , (23) 

V ZMI\ ev t s 

independent of the order 2k, and much larger than the estimate (Ell) 



We also discuss in Appendix D2 the non-Gaussian character of the fluctuations of the estimates v n {2k}, due to 
their non-linear relations with the cumulants, Eqs. (fT^). In particular, statistical fluctuations may result in a "wrong" 
sign for the cumulant, in the sense that v n {2k} 2k defined by Eqs. @ is negative. If this happens, the flow clearly 
cannot be estimated from the corresponding cumulant. 



E. Non-unit weights 



In the definition of the generating function, Eq. (|^), each particle was given the same weight. A more general form 



is 



M 



3 = 1 



(24) 



The weight Wj can be any arbitrary function of the rapidity y of the particle, its transverse momentum px, its mass. 

In order to obtain the highest accuracy on the flow measurement, Wj should be chosen proportional to the flow 
itself, as shown in E?|] (see also pTj| , pT[ ): the ideal weight is w(px,y) oc v' n (j>T,y), which is intuitively clear: one must 
give higher weights to particles with stronger flow. This is also the best choice if one uses the standard method, 
involving the determination of the reaction plane. 

Ideally, the flow analysis should be performed twice: a first measurement of the flow can be done with reasonable 



guesses for the weights; measuring differential flow as a function of y and px for various particles (see Sec. ill), the 
values obtained can be used as the new weights in a second, more accurate analysis. This is the procedure recently 
followed in |p6| . 

A variety of weights have been used in analyses of the directed flow V\ . Since v\ changes sign at midrapidity, the 
weight must be an odd function of the rapidity in the center-of-mass frame. Most often, the weight is simply the 
sign of y, with jl7) or without j52| a gap at midrapidity. A linear dependence in y was used in ]53|-|55|]. The latter 
choice is better, since v\ is itself linear near midrapidity. The transverse momentum dependence of w is most often 
linear, as in the original paper p7| . Unit weights, independent of pr, are also widely used [^1 28 5(||57|]. They are 
convenient, because no particle identification is required. However, since v\ is linear in pt, (at least at low pr plfl) , 
the original choice w cx pt is likely to give more accurate results, although the opposite conclusion was reached in jg6|. 
At intermediate energies, one can in addition choose a weight proportional to the mass of the particle, to take into 
account the fact that nuclear fragments flow more than protons [ |58|]59| ] . 

Unlike directed flow, elliptic flow is an even function of the center-of-mass rapidity: therefore the weights are 
usually chosen independent of rapidity. The weights are either independent of transverse momentum P Jl^ , |28| , p7| , p8[ 
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or proportional to p\ i MH^]^ ]- The latter choice is more appropriate if px is measured, since v 2 is proportional to 
At ultrarelativistic energies, however, v 2 is almost linear in pt above 100 MeV @J6l|. A better 

Po, with po ~ 100 MeV: this weight, quadratic at low pt and linear 



p\ at low |51 



choice may be for instance w(j>t) — yVr 
at high px, reproduces roughly the measured pt dependence of v 2 . 

With a non-unit weight w, the following modifications should be applied: 



• In the relations between the cumulants and the flow, i.e., Eqs. (|17]) for integrated flow, Eqs. (|34|) and fl3q ) 
for differential flow, and the equations of Appendix O, v n now stands for (we m ^"*«)). Comp aring with the 
previous definition, Eq. (pj) , the flow is now weighted by w, as can be expected from the definition of the 
generating function, Eq. (|24|). 

• In the formulas giving the statistical errors, derived in Appendix O, on e must in addition replace M by M / (w 2 "}- 



As a consequence, the resolution x appearing in Eqs. (D21), (D22), (D35) and (D36) now stands for x 



M (wv n ) / (w 2 ). 

• In the interpolation formulas of Appendix 



w, for instance by \J (w 2 ). 



the value of r [Eq. (Bf )] should be scaled by the typical value of 



III. DIFFERENTIAL FLOW 



When one has measured the flow integrated over phase space, the next step is to move on to the "differential 
flow" analysis, i.e., the measurement of flow in a narrower phase space window. We call a particle belonging to the 
window of interest a "proton" (although it can be anything else). We denote by ?/> its azimuthal angle, and v' n its flow 
harmonics (the so-called "differential flow"), v' n = ^e 1 "^ - *")). The particles used to estimate the integrated flow v n 
are named "pions" . 

In order to measure the differential flow of the protons, we correlate their azimuth ip with the pion azimuths 
4>j. Once the integrated flow v n is known, one can reconstruct the differential flow v' n from this correlation, and 



also higher harmonics v' 2rr 



u 3m 



etc. It follows that differential elliptic flow, v 2 , can be reconstructed using either 



integrated directed flow or integrated elliptic flow. Following |lq| , we denote by v' . the differential flow v' p measured 
with respect to integrated flow v n , where p is a multiple of n. At intermediate energies, one usually measures v' 2 



2I| while Vy 2 is more accurate at ultrarelativistic energies where v\ becomes very small |28 



2/1 



The differential flow is reconstructed b y tak ing the cumulants of azimuthal correlations between the proton and 
the pions. These are constructed in Sec. [II A by means of an appropriate generating functi on. T he subtr action of 
"autocorrelations", in the case when the proton is one of the pions, is briefly discussed in Sec. [II B. In Sec. [II C, we 
derive the relations between the cumulants and the differential flow, 



As in the case of integrated flow, cumulants 
of different orders yield different estimates of v' p / n ■ The optimal choice is that which minimizes uncertainties, discussed 



"p/r. 



in Sec. HID 



A. Cumulants 



To measure a proton differential flow v' p , we first construct a generating function of measured azimuthal correlations 
between the proton and pions. This function is the average value over all protons of e tp ^G n (z), where G n (z) is the 
generating function defined by Eq. (|J), evaluated for the event where the proton belongs. Note that the average 
procedure is not exactly the same as when studying integrated flow. One must first average over the protons in a 
same event [i.e., with the same G n (z)]; then, average over the events where there are "protons" (if the "proton" is 
some rare particle, or if the phase space window is small, there may be events without a proton). 

Expanding in power series of z and z* , one obtains 



(e ilnp G n (z)) = (e 1 ^) + z (V^-™* 1 ^ + z* / e ^ +n<?il ^ + 



(25) 



This generates measured azimuthal correlations between the proton and an arbitrary number of pions. The generating 
function of the cumulants is a complex- valued function of the complex variable z: 



Vp/niz) 



= (e^G n (z)) 

~ (Gn(z)) ■ 



(26) 
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where (G(z)) denotes an average over all events, as in Sec. [II B . The cumulants are by definition the coefficients in 
the power-series expansion of this function: 



V p/n {z) = J2 ~j^J[~ ^e ip * +iB {* 1+ - + **-**+ 1 — ^ . (27) 

k,l 

The physical significance of these cumulants is the same as for the cumulants used in the analysis of integrated flow. 
They eliminate detector effects and lower-order correlations, so that only the direct correlation between k + I pions 
and the proton, of order M~ k ~ l , and the correlation due to flow remain. If the proton is not correlated with the pions, 
for instance, Eq. (26) gives V p / n (z) = (e* p ^), independent of z, and all cumulants with k + I > 1 vanish according to 
Eq. 



27). In the more general case when there are correlations, expanding Eq. fl26p to order z and identifying with 
Eq. (27), one obtains 

| e *W-^i)| = / e i(p<l>-n<f>i)^ _ (giP^ ^e - " 1 ^ 1 ) . (28) 

This cumulant is analogous to Eq. (||), and can be interpreted in the same way. 

If the detectors used to measure protons and pions are perfect, simplifications occur: first, the cumulants defined in 
Eq. (^) are real. To show this, we use the property that if the detector is perfect, the probability that an event occur 
is unchanged when one reverses the sign of all azimuthal angles (i.e., tf> — > —?p and <pj — ► —<f>j), therefore T> p / n (z) is 
unchanged under this transformation. Now, the transformation <pj — > —<f>j amounts to changing z into z* in G n (z), 
according to Eq. @. Thus, Eq. (|||) can also be written V p / n (z) = {er ip ^ G n (z*)) / (G n (z*)). Comparing with the 



original definition, Eq. (26), z has been changed to z* and iji to —tp. Now, changing ip to —ip in Eq. (|26|) amounts 
to taking the complex conjugate of D p / n (z), since G n (z) is real. One finally obtains T> p / n {z) = T>*^ n (z*), from which 
one easily deduces that the coefficients in the expansion (|27l ) are real. They are in general complex with a realistic 
detector, but only the real part is relevant. 

Furthermore, most cumulants vanish if the detector is perfect. In order to see this property, we shift all angles by 
the same quantity a, which does not change the probability of the event. The a ngles of the pions are changed to 
4>j — a, which amounts to changing z into ze ma in G n {z), as explained in Sec. |ll B| . Similarly, the angle of the proton 
ip is changed to ip — ► ip — a, so that e lp ^G n (z) becomes e ip v e~ tpa G n (ze ina ) . Averaging over a gives unless p is 
a multiple of n, which is the case we are interested in. Writing p — mn, the only terms which remain in the power 
series expansion of (e lp ^G n (z)) are the terms in z * k z k + n \ To obtain the generating fun ction of cumulants T> p / n (z), 



one must divide by (G n (z)). Since this quantity depends only on \z\, as explained in Sec. [IB, the only nonvanishing 
cumulants in Eq. ( p^ ) are those with I = k + m. Other cumulants are physically irrelevant. 
Finally, the relevant quantities are 



{2fc + m + l} = 5R // e *»("nM-*i+-+**-**+i— -*■*+») \\ . ( 2 g) 



where 5ft denotes the real part, and the notation {2k + m + 1} means that the cumulant involves correlations between 
2k + m + 1 particles (one proton and 2k + m pions). 

When there is no flow, d mn/n {2k + m + 1} is of o rder M~ 2k ~ m . Flow gives a contribution to this cumulant, 
proportional to v' p v^ +m , which is calculated in Sec. Ill Cj . If this is the dominant contribution, one obtains an 



estimate of the differential flow v' p from the cumulant d p / n {2k + m + I}, using a previously determined value of the 
integrated flow v n . This estimate will be denoted by v' p ^ n {2k + m + I}. 

Analytical expressions of higher order cumulants, deduced from Eq. (E7j), are lengthy. As in the case of integrated 
flow, the simplest way to extract them consists in tabulating the generating function T> mn / n (z), and then computing 
numerically the coefficients of its power-series expansion, through interpolation formulas which are given in Appendix 



C2 



B. Autocorrelations 



When studying the azimuthal correlations between the "proton" and the "pions", the "proton" must not be one 
of the "pions" , otherwise trivial autocorrelations would appear. This problem is well known in the standard flow 
analysis p7| : in order to avoid it, one excludes the particle under study (the "proton") from the definition of the flow 
vector used to estimate the reaction plane, which usually involves all the other particles (the "pions"). Here, if the 
proton is one of the pions, one simply removes its contribution by dividing G n (z) by 1 + (z*e m ^ + ze~ m ^)/M in the 
numerator of Eq. ( p6| ) . The generalization of this procedure to non-unit weights is straightforward. 
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C. Contribution of flow to the cumulants 



Let us now calculate the contribution of flow to the cumulants d mn / n {2k + m + 1}. As in Sec. II C , we neglect 
nonflow correlations and assume a perfect detector for simplicity. Under these assumptions, we can compute the 
generating function of the cumulants T> p i n (z). We first average over all events with a fixed orientation of the reaction 
plane $_r, and over the protons in each single event: 

(e**!**) =v' p e ip *». (30) 



Hence, Eq. (|26| ) becomes 

Vp/niz) = — (G~Tzj) Vp ' ' 

where the denominator is given by Eq. (jl^), and (G n (z)\^n) by Eq. (|l3|). The numerator vanishes unless p is a 
multiple of n, i.e., p — mn with m integer. Integrating over one then obtains: 



z *k gk+m 



Z 



~I m (2v n \z\) l^r J , (32) 

where, in the last identity, we have assumed that M is large, so that Ml / (M — in — 2k) ! ~ M 2k + m and we may extend 
the sum over k to infinity. Equation ( |3l| ) gives 

_ /™(2|£K) / z \ m , 

^mn/n{Z) - j , v 111) V' 



Jo(2|z|«„) Vkl 

This equation can be expanded in power series of z and z* . The coefficients of the power-series expansion are real, 
as expected from the discussion in Sec. Ill A. Comparing with Eq. (|27|), cumulants with I ^ k + m vanish, which was 



also expected; cumulants with I — k + m, which are the d mn / n {2k + m + 1} introduced in Eq. (29), are proportional 
to v' mn . We thus obtain estimates of the differential flow v' mn , which we denote by v' mn / n {2k + m + 1}, from the 
measured cumulants. For m — 1 (the proton is correlated with pions in the same flow harmonic) these estimates are 
given by 

</„{2} = d n/n {2}/v n (34a) 
v' n/n {4} -d n/n {4}/vl (34b) 

while for m = 2 (useful when measuring differential elliptic flow v' 2 with respect to the integrated directed flow v\ ) 

w 2n/«{3} = d 2n /n{3}/vl (35a) 
«an/»{5} = -d2n/n{5}/(2<). (35b) 

As in the case of integrated flow, a non-perfect acceptance may induce interference between the various harmonics 
v' n , v' 2n , etc., modify ing Eqs. ( |34| ) and (p5|). The corresponding relations between the cumulants and the flow are 
derived in Appendix C2. We have taken into account the possibility that integrated and differential flows may be 
measured using detectors with different azimuthal coverages, which is often the case in practice (see for instance pl[). 
In particular, if the detector used for integrated flow is perfect, we show that no correction is required for v' n , whatever 
the detector used for differential flow may be, which is intuitively obvious: the only difference when using a smaller 
detector for the reconstruction of differential flow is then a loss in "proton" multiplicity, resulting in higher statistical 
errors, which we now discuss. 



D. Errors 



As in Sec. II D, we now evaluate the contributions of nonflow correlations and statistical fluctuations to the cumu- 



lants. This will allow us to determine the optimal cumulant order to be used, which minimizes the total uncertainty 
on v'. 
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As discussed in Sec. [II A, nonflow correlations give an unknown contribution of order M~ 2k ~' m to the cumulant 
d p /„ {2k + m + 1|, w hich must be compared with the contribution of flow, proportional to v' mn {2k + m + 1} v„ +m as 
shown in Sec. Ill C . The systematic error on v' mn , thus reads 



<m/«{2* + m + 1} - v' mn/n = O {(Mv n )- 2k - m ) ■ (36) 
According to Eq. (EG), this systematic error becomes smaller and smaller as k increases: the same behavior was 



observed for the systematic error on the integrated flow in Sec. II D . 

The order of magnitude of statistical errors can be estimated in the same way as for integrated flow. The cumulant 
d p / n {2k + m + 1} involves correlations of the proton with 2k + m pions belonging to the same event. There are roughly 
M 2fc+m ways (if M is large enough) to choose 2k + m pions among M. Denoting by N 1 the total number of protons 
involved in the analysis, there is a total of M 2k+m N' subsets of particles involved in evaluating the cumulants. The 
resulting statistical error is therefore 

v' /n {2k + m + 1} - v' /n = O ( . 1 = ) }, . (37) 

As was the case for integrated flow, the statistical error usually increases with the order of the cumulant, while the 
systematic error decreases. Therefore, the cumulant order which gives the best compromise is that for which statistical 
and systematic errors are both of the same order of magnitude. Equating the right-hand sides of Eqs. ( |36| ) and (|3^), 
one obtains the optimal cumulant order |47f| : 

In N' 

(2fc + m + l) opt ^l + — . (38) 



Statistical errors are evaluated in detail in Appendix p 3| . We show that the simple estimate (|37j) is correct only 
if x = v n \/~M <C 1. In this limit, different estimates v' mn ^ n {2k + m + 1} and v' mn ^ n {2l + m + 1} with k ^ I are 
uncorrelated. When ^> 1, the correlation becomes strong and the standard error is approximately 

*«,«/«{2* + m + l}) stat - -4=, (39) 



independent of the order 2k + m + 1 , and much larger than the estimate (J37]) . 



IV. COMPARISON WITH OTHER METHODS 



The method proposed in this paper is closely related to our first cumulant-based m ethod p7| , since it is aimed at 
correcting some of the latter's limitations. This is discussed in Sec. 1VA. Then, in Sec. [VB, we compare our method 
with the two-particle correlation technique |3(| and with the subevent method |l^,|l8| ■ 



A. Comparison with our previous method 

The cumulant expansion proposed in Efl] w as based on the flow vector rather than on particles themselves. The 
flow vector is defined for each event by Jl7|]l^] 



1 M 

— V e m ^ . (40) 



3 = 1 

A generating function was then defined: 

g (z)= (e zQ * +z * Q «). (41) 

The expansion of this generating function generates all the moments of the distribution of Q, that is, (Q k Q* 1 )- One 
easily shows that Oo(z) is closely related to the generating function used in the present paper. Using the identity 
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valid for large M, we may rewrite the generating function (|^) as 



(M M 
j=i i=i 

= exp ( 4=Q„ + -^=Q* n ) . (43) 



M VM 



Thus the average over events (G n [z)) coincides with Qo(z/yM). This shows that both methods are equivalent in the 
large M limit, up to a rescaling of the variable z. 

The difference between the two approaches can be easily understood by expanding the two generating functions in 
powers of z and z*. To second order, for instance, Eqs. ( |40| ) and ( ff"H ) give 



••• ' 77(2. e^"**J ) + — ( ^ e -^+«*J ) + — ( ^ e m ^+<^ ) + • • • (44) 




while our new generating function gives (see Eq.(^)) 

{Gn{z)) = • • • + g e^"*») \ + ^ (V e -(^)\ + ^ e^^A + • • • (45) 

The essential difference is the restriction j 7^ k in the sums: our new method is free from the autocorrelations 
corresponding to the terms with j = k in Eq. (Q). This remains true to higher orders in z and z*. 

Autocorrelations have two effects: in the first term of Eq. (fhH), they give a constant, trivial contribution which 
must then be removed. This fixes the choice of the weight 1/V M in the definition of the flow vector, Eq. (pli|). With 
another weight, the contribution of autocorrelations would depend on M, and it would not be easy to subtract them 
when M is allowed to fluctuate. With the new generating function, we are free to choose another weight, and we 
show in Appendix [A that the weight 1/M [Eq. (||)] gives more accurate results when M is allowed to fluctuate. 

In the second and third terms of Eq. (p4|), autocorrelations create terms in e ±2m * J which interfere with higher 



harmonics. This was the main limitation of the method exposed in 47 1. Eliminating all autocorrelations represents 
a major improvement. In particular, our method should enable to measure directed flow at RHIC, if any. 

Finally, let us comment on our definition of the cumulants through the generating function C n (z) in Eq. (Q). This 
definition ensures that cumulants of order 2 and higher vanish if particles are uncorrelated [see the discussion following 
Eq. (g)]. In the limit when M is large, one recovers C n (z) = In (G n (z)), in agreement with the standard definition of 
cumulants in probability theory flO, and with the definition we adopted in W7|. 



B. Comparison with standard methods 

In [p0|| , it was proposed to analyse flow using two-particle azimuthal correlations. More specifically, one defines 
by P C orr(A</>) [resp. -P U ncor(A^>)] the distribution of the relative angle A<j) = <j>\ — 4>2, where 4>i an d 4>2 are any two 
particles belonging to the same event (resp. belonging to two different events). One then constructs the ratio 

C(Afl= p Pc ° rr( f A ^ v (46) 

This so-called "mixed event" technique enables the extraction of the physical correlations between the particles, 
eliminating the effects of an uneven detector acceptance. Neglecting nonflow correlations, one has in general 

+00 

C(A0)= v 2 J nA t (47) 

n— — 00 

so that the Fourier expansion of the measured correlation function C(A<j)) simply yields the integrated flow v n . 
Similarly, if one replaces <j>i with the azimuthal angle ip of a particle in a narrow phase space window, v\ in Eq. (^7j) 
is replaced with v' n v ni where v' n is the differential flow of the particle of interest. 
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The values of v n and v' n obtained with this method coincide with the values we obtain from the cumulant of order 
2, which are denoted by v n {2} and v' n {2} in this paper; while we do not need mixed events, we must apply correction 
factors to our re cons truc ted ^n {2} and ^{2} if the detector does not have full azimuthal coverage, as explained at 



the end of Sees. II C and QIC. With the mixed-event technique, this correction is not required |33| . 

The essential limitation of the two-particle correlation technique is that it is not possible to separate flow and 
nonflow correlations: therefore the results may be strongly biased by nonflow correlations, which arc eliminated in 
our method by the construction of higher order cumulants. 

The same limitation applies to the much more widely used subevent method p7|| . This methods involves two steps: 
in each event, one constructs the flow vector ( flO| ) to estimate the orientation of the reaction plane <&r, and a study 
of the azimuthal correlation between the flow vectors of two randomly chosen "subevents" yields the accuracy of this 
estimate (the so-called "reaction plane resolution"). This first step amounts to measuring the integrated flow. Then, 
in order to measure differential flow, one studies the azimuthal correlation between a single particle and the flow 
vector. Since the latter involves a summation over many particles, the correlation between a single particle and the 
flow vector is much stronger than the correlation between two single particles, which is probably the reason why this 
method is used so often. However, the relative weights of nonflow and flow correlations is the same [^J as in the 
much simpler two-particle correlation technique discussed above, so that both suffer from the same limitations. 

Flow and nonflow correlations could in fact be distinguished in the subevent method, through a more detailed study 
of the azimuthal correlation between subevents. Indeed, flow and nonflow correlations yield a different shape of the 
distribution of the relative azimuthal angle A<!> between subevents, as shown in |63|]. If correlations are only due to 
flow, then the distribution of A$ is universal, and depends only on the resolution \ = v n \f~M which characterizes 
the reaction plane resolution fli^| . A comparison of the calculated distribution with experimental data was recently 
done at energies of 250 MeV per nucleon H. The agreement is perfect, which shows that the observed correlations 
are dominated by flow at these energies. To our knowledge, no such comparison has been carried out so far at 
ultrarelativistic energies. 

V. RESULTS OF MONTE-CARLO SIMULATIONS 

We have performed various Monte-Carlo simulations to check the validity of the procedures explained in this paper. 
In each simulation, A^ cv t s events are simulated; in a given event, the reaction plane is chosen randomly, then 
"pions" and "protons" are generated according to the respective distributions 

dN 

— oc 1 + 2vi cos(> - $r) + 2v 2 cos (2(0 - $#)) (48) 



and 

dN' 



, , cx 1 + 2v[ coslip - $r) + 2v' 2 cos (2(ip - $ R )) (49) 
dip 



respectively. We then reconstruct the integrated (Sec. V A) and differential (Sec. VB) flows following the procedures 
exposed in Sees. || and III. 



A. Integrated flow 

In the first set of simulations, we generated N evts = 10 5 events with M = 200 pions emitted with an integrated 
elliptic flow V2 = 6%, which we then tried to reconstruct. No integrated directed flow was simulated. We first assumed 
a perfect detector. The optimal cumulant order, defined by Eq. (|22j), is 2fc opt = 4.2. The estimate t>2{4} derived from 
the fourth-order cumulant is thus likely to give the best compromise between systematic and statistical errors, but 
we also calculated the estimates v 2 {2} and «2{6}. 

The results are presented in Table |. The error bars are statistical only. They are asymmetric for «2{6}: this reflects 
the fact that the statistical fluctuations of v n {2k} are not Gaussian when the error is large, as explained in Appendix 



D2. When working with a fixed multiplicity M = 200, the reconstructed values coincide with the theoretical value 
within expected statistical errors. Note that the statistical error on V2{6} is only slightly larger than that on t>2{4} in 
this case. When the multiplicity M is randomly chosen between 150 and 250, the reconstructed flow deviates from the 
theoretical value by more than two standard deviations, but the accuracy is still good, and moreover the statistical 
errors were calculated assuming a fixed detected multiplicity, while the real errors are probably larger. 
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v 2 {2} 


« 2 {4} 


w 2 {6} 


full acceptance, M = 200 


6.01 ±0.02 


6.03 ± 0.04 


6.08 ±0.05 


full acceptance, 150 < M < 250 


6.01 ±0.02 


6.11 ±0.04 


6.13 ±0.05 


nonflow correlations, M = 200 


9.21 ±0.02 


6.16 ±0.04 


6.19 ±0.05 


"bow tie" acceptance, M — 100 


5.99 ± 0.04 


6.09 ±0.12 


R n A +0.19 

6.34 _ Q 22 


"bow tie" acceptance, 75 < M < 125 


5.85 ±0.04 


5.40 ±0.12 


5 74 +0.19 
d -' -0.22 



TABLE I. Reconstruction of integrated elliptic flow. iV cv ts = 10 5 events with vi = and v 2 = 6% were generated in each 
simulation. The three columns give the values of v?, in %, reconstructed using cumulants to order 2, 4, and 6, respectively. 



In order to check the ability of our method to eliminate two-particle nonflow correlations, the latter were simulated 
by emitting particles in pairs, where both particles in a pair have the same azimuthal angle. 100 pairs were emitted in 
each event, resulting in a multiplicity AI = 200. The reconstructed value with the cumulant of order 2 , i.e., without 
removing two-particle nonflow correlations, is 50% too large; the standard methods would give similar results. On 
the other hand, the values obtained using higher order cumulants are much closer to the theoretical value; they are 
beyond statistical error bars, but this is not surprising since error bars were calculated assuming M = 200 independent 
particles, while the effective multiplicity here is rather 100, which results in larger fluctuations. 

In summary, the results obtained so far show that v 2 {4} is to be preferred here: the statistical error on ^{6} is 
(slightly) larger, while nonflow correlations may give a large, uncontrolled, contribution to i> 2 {2}. 




FIG. 1. Schematic picture of the PHENIX detector at RHIC. The shaded area indicates the azimuthal coverage of the 
detector, while the darker area at the right of the figure corresponds to a small time-of-flight detector with an extension of 45°, 
which can be used to measure the differential flow of identified particles. 



To test the validity of our acceptance corrections, we then did simulations with the "bow-tie" detector schematically 
represented in Fig. El which mimics the azimuthal acceptance of the PHENIX detector [Q: particles are detected 
only within two quadrants of 90 degrees each, with 100% efficiency. Since only half the particles are detected with this 
detector, the value of M chosen in this simulation was half the value chosen above for a perfect detector. The optimal 
cumulant order is now 2fc opt = 4.5, so that the fourth-order cumulant is still to be preferred. Since the azimuthal 
coverage is only partial, Eqs. (|l7|) relating the cumulant to the flow are no longer valid. In the case of elliptic flow, 



they are replaced by Eqs. (C7). These formulas involve the Fourier coefficients a p of the acceptance function, defined 
by Eq. (|l8|). With the acceptance depicted in Fig. [j], a simple calculation gives a 2 — 2/n and a\ = a% = = 0. 



Inserting these values into Eqs. (C7), one finds 

c 2 {2}~ 0.518 v 2 {2} 2 , 

c 2 {4}~ -0.384 v 2 {4} 4 . (50) 

An interesting feature of this bow-tie acceptance is that there is no interference between v\ and v 2 > because all a p with 
odd p vanish. However, this is not the case in general (see Appendix ^J). The results are given in Table || Note that 
the estimate of statistical errors was done assuming a full acceptance, so they may be underestimated here. When the 
multiplicity M is fixed, the reconstructed values agree with the theoretical values within eror bars. If the acceptance 
correction had not been applied, the reconstructed values would be below 5%. 

If M is allowed to fluctuate, the discrepancy between reconstructed and theoretical values is much larger than 
statistical errors, and also much larger than in the case of a full acceptance. As explained in Appendix |A], a fluctuation 
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SM in multiplicity induces some errors, which scale like SM 2 and are comparatively larger when the acceptance is not 
good. In order to avoid this effect, the analysis can be don e by selecting randomly a subset of the detected particles, 
with fixed multiplicity, as explained at the end of Sec. |IIB|. 





V 2 = 


v 2 = 6 % 


full acceptance, M = 200 


2Q1 +0.21 
Z ' U1 -0.32 


2 19 +0 - 21 

-0.32 


"bow-tie" acceptance, M = 200 


, 07 +0.21 
l" 87 -0.32 


1 00+O.2I 



TABLE II. Reconstruction of a theoretical directed flow vi —2%. N ev ts = 5 x 10 5 events were generated in each simulation, 
with M — 200 detected particles in each event. The table gives the values of i>i{4} in %. 



We then performed a second set of simulations, with both directed and elliptic flow, in order to test possible 
interferences between the two. Such uncontrolled inter ferenc es were the main limitation of our previous cumulant 
method Gt), but are avoided here, as explained in Sec. IV A. Inteference of a d ifferent kind may in fact still occur 



when the detector has partial azimuthal coverage, as explained in Appendix C 1 , but they are under control. Events 
were generated with a small directed flow v% — 2%, and with an elliptic flow V2 = or V2 = 6%. The goal was to 
reconstruct v\. Since 2% is a very small value, a larger number of events were generated than in the previous set of 
simulations. With the values of iV cv ts and M given in Table 0, the optimal cumulant order is still 4: all the results 
presented in Table [n| are reconstructed from t he fo urth-order cumulant. In the c ase of the bow-tie acceptance, the 
relation between «i{4} and the cumulant, Eq. (17b), must be replaced with (C6b), which gives 



ci{4} ~ -2.785 vi{4}" 



(51) 



The reconstructed value of v\ is always compatible with the theoretical value v± — 2% within the statistical error, even 
with the bow-tie acceptance. This is a major improvement on our previous method p7|, where the reconstruction of 
Vi would have been impossible if V2 = 6%, even in the case of a perfect detector. 



B. Differential flow 



We then performed various simulations to test our reconstruction of the differential flow. In a first set of simulations 
(see Table III), only elliptic flow was generated. Events were generated with v' 2 = 1)2 = 6%, and we tried to reconstruct 
v' 2 - As explained in Sec. p]| , the first step in the analysis is the reconstruction of V2- Since no integrated directed 
flow vi is present, v 2 can then be reconstructed only with respect to V2: this is denoted by Vy 2 in Sec. IIIJ One 



can reconstruct it from the lowest order cumulant. This yields the estimate Vy 2 {2}, equivalent to the standard flow 
analysis. One can also reconstruct the higher order cumulant Vy 2 {A}^ which eliminates nonflow correlations between 
the proton and the pions. 



pions 


protons 


«2/2{2} 


«2/ 2 {4} 


full acceptance, M = 200 


full acceptance 


5.96 ±0.15 


6.10 ±0.44 


"bow-tie" acceptance, M = 100 


"bow-tie" acceptance 


6.07 ±0.20 


5.83 ±0.99 


"bow tie" acceptance, M = 100 


45° acceptance 


6.02 ± 0.20 


5.71 ±0.99 



TABLE III. Reconstruction of differential elliptic flow with respect to integrated elliptic flow, i.e., u 2 /2- ^ n eac h simulation, 
TV' = 5 x 10 5 protons were generated with v' 2 = 6% and v[ — 0. The integrated flow was also U2 = 6% and v\ — 0. The last 
two columns give the values of v' 2 (in %) reconstructed from the cumulants to order 2 and 4, respectively. 



In the first simulation, we assumed that the detector was perfect, both for integrated and for differential flow. With 
the values of the parameters given in Table III , the optimal cumulant order, defined by Eq. (|38|), is (2fc±m±l) op t = 3.5, 
i.e., i>2/2{4} is to be preferred. Nevertheless, we also include the lowest order estimates «2/2{2}. Since no nonflow 
correlations are generated between the proton and the pions, this estimate should be good. Indeed, both estimates 
are found to be in very good agreement with the theoretical value. 
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We then turned to the case whe n the detector no longer has a perfect acceptance. In this case, corrections must be 
appl ied, as explained in Appendix C2. For the integrated flow, we used the "bow-tie" detector of Fig. |l|. As in Sec. 



V A. the multiplicity of each event with this detector is M — 100, instead of M = 200 with the full acceptance. This 
modification does not affect the optimal cumulant order, which is still 4. For differential flow, two different detectors 
were considered. First, the same bow-tie acceptance as for pions, i.e., the relevant Fou rier c oefficients arc a 2 = 2/ir 
and a[ = a' 3 = a' 4 = 0. In this case, the relations between the cumulants and flow, Eqs. (|C15| ), are 

d2 /2 {2}~ 0.595 v 2 v' 2/2 {2}, 

d 2 /2H}- -0.406 v*v' 2/2 {4}. (52) 

Secondly, we assumed the protons were detected with a small 45° detector (see Fig. [t]) . The Fourier coefficients of its 
acceptance function A'(ip) are a' = 1, a' p>1 = 8sin(p7r/8)/(p7r), which yields the relations 

d 2/2 {2}~ 0.337 v 2 v' 2/2 {2}, 

d 2/2 {4}~ -0.181 v 3 2 v' 2/2 {4}. (53) 
Note that the last correction is more than a factor of 5, compared to a full acceptance. In all cases, the results (Table 



III) are in excellent agreement with the theoretical value. This shows that with our method, it is possible to measure 



differential flow with any detector provided the relevant acceptance corrections are performed. 



pions 


protons 


«2/l{3} 




« 2/2 {4} 


full acceptance, M = 200 


fall acceptance 


6.47 ±0.32 


5.81 ±0.15 


6.03 ± 0.44 


"bow-tie" acceptance, M = 100 


"bow-tie" acceptance 


5.08 ±0.54 


5.80 ±0.20 


5.51 ±0.99 


"bow tie" acceptance, M = 100 


45° acceptance 


5.20 ±0.54 


6.81 ±0.20 


7.95 ±0.99 



TABLE IV. Same as Table En, except that a directed flow vi = v[ = 6% was generated in addition to the elliptic flow. The 
third column gives the differential elliptic flow recontructed with respect to integrated directed flow, to lowest order. 



We then did a second set of simulations. The only difference with the first set is that directed flow is also generated, 
with the same magnitude as elliptic flow, i.e., v[ = v\ = 6%. In this case, it is also possible to measure differential 
elliptic flow with respect to integrated directed flow, v 2 /r We recontructed only the lowest order estimate t;^^}, 
which is also given by the standard flow analysis. Note that it is already insensitive to two-parti cle no nflow correlations. 
In the case when the bo w-tie acceptance applies to both integrated and differential flow, Eq. (C14) does not give any 
corre ction , so that Eq. (35a) still holds. When protons are measured with a small 45° detector, on the other hand, 
Eq. (IClJ) gives 



{3} ~ 1.258 vlv' 2/1 {3} 



(54) 



Results are given in Table [V. With a full acceptance, the reconstruction of both v' 2 , x and v' 2 / 2 is good, within the 
expected error bars. With a partial acceptance, however, discrepancies are much larger. In particular, the value of 

' „rUY, „ AKO „„„„„ 
2/2 



v' n /n with a 45° acceptance differs by four standard deviations from the theoretical value 



VI. CONCLUSIONS 



We have presented a new, general method for analyzing the flow. This method can be used by all heavy-ion 
experiments which implement the standard analysis of Danielewicz and Odyniec ]l7|] , since it uses the same input, 
namely the azimuthal angles of outgoing particles. As in the standard analysis, one proceeds in two steps. One first 
reconstructs the average value of the flow over phase space, with an appropriate weight: this "integrated flow" v n is 
related to the reaction plane resolution in the standard analysis. One may then perform detailed analyses in narrower 
phase-space windows, i.e., measure "differential flow" v' n . The strong points of our method are the following: 

• It systematically eliminates azimuthal correlations which are not due to flow. 

• One obtains several different estimates of v n and v' n ] comparison between these estimates provides a useful 
consistency check. 
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• The detectors need not have full azimuthal coverage. 

Let us comment on these three points in more detail. 

Nonflow correlations are eliminated by means of a cumulant expansion of multiparticle azimuthal correlations. This 
cumulant expansion applies to both integrated and differential flow. In the case of integrated flow, one constructs a 
series of cumulants c n {2k}, where k is an arbitrary positive integer. c n {2k} involves correlations between 2k particles, 
but is insensitive to correlations involving less than 2k particles. Now, flow is by definition a collective phenomenon, 
which induces correlations between all the produced particles: as k increases, its relative contribution to c n {2k} also 
increases. 

Using the cumulant c n {2k}, one obtains an estimate of v n , denoted by v n {2k} in this paper. The lowest order v n {2} 
corresponds to the standard flow analysis: it coincides with v n if all correlations are due to flow. The nice feature 
of our method is that higher order cumulants provide us with other, independent estimates of v n , which can be used 
as a consistency check. The reader might believe that constructing cumulants of six-particle correlations requires a 
huge statistics, and is practically impossible. This is not true. Statistical errors on higher order estimates v n {2k} 
are discussed thoroughly in Appendix [d|. They strongly depend on the strength of the flow itself, and more precisely 
on the parameter x = v n VM, which characterizes the resolution of the reaction plane reconstruction, A^, in the 
standard analysis. If \ = 1) corresponding to a resolution = 45°, which is achieved in many experiments, the 

statistical error on both w„{4} and u„{6} is only 60% higher than the error on v n {2}, as can be seen in Fig. ||. 

The same discussion can be repeated for differential flow. In the case of differential elliptic flow v' 2 , however, the 
standard analysis also provides us with two different estimates of this quantity, depending on whether the reaction 
plane is measured using directed {v' 2 /i) or elliptic (vy 2 ) how. The first of these two estimates, Vy^, which involves 
correlations between three particles, is insensitive to two-particle nonflow correlations, even in the standard analysis. 
In this paper, we obtain several estimates for both v 2 ^ and v' 2 ^ 2 , which can be used as further consistency checks. 

Finally, we have also included a detailed discussion of acceptance corrections, which allow us to work even if the 
detector has only partial azimuthal coverage. One may argue that the corresponding formulas, given in Appendix 
, are very heavy. However, these formulas are required only if the acceptance is far from isotropic: corrections 
corresponding to weak inhomogeneities of the acceptance are automatically taken care of by the cumulant expansion 
itself. Furthermore, the correction factors need only be calculated once for a given detector, and it is a simple 
calculation. When they are applied, accurate results can be obtained even with a very poor acceptance, as illustrated 
by the simulations performed in Sec. JM. Although several methods have been proposed to correct for detector 
inefficiencies in the standard analysis JL8|, we do not know of any systematic study of their limits. 
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APPENDIX A: CUMULANTS OF MULTIPARTICLE AZIMUTHAL CORRELATIONS 



In this Appendix, we justify the definition of the cumulants from the generating function (Q), by discussing a few 
examples. 

We begin with the second order cumulant c n {2} = ie m ^ >1 ~^\. Inserting Eq. (JsJ) in expression ^ and expanding 
to order zz*, we obtain from Eq. (pi) 



^ e M*i-*0^ = ^L-i ((eM*x-h)} - (e'"*)(r'"*}) . (Al) 

If particles are uncorrelated, i.e. (V"^ 1- ^ 2 * 1 ) = (e™^ 1 ) (e~ m< ^ 2 ). the cumulant vanishes, as it should. In the limit 
when M ^> 1, one recovers definition (||). 

The standard definition of the cumulan ts i s C n {z) = ln(G„(z)), rather than Eq. (Q) |52j. However, with this 



definition we would obtain, instead of Eq. ( Al), 

M - 1 



e m(0i-0 2 )\\ _ 



M 
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M-l 

M 



(A2) 



If particles are uncorrelated, the last term in the equation remains if the acceptance is not perfect. This may result 
in large e rrors when the flow is small. 

In Sec. II B , we also considered the possibility that the multiplicity M be not strictly the same for all events. If 
it fluctuates around an average value (M), then the generating function (^) m ust a lso be averaged over M, and M 
must be replaced by (M) in the definition of the cumulants, Eq. (M). Equation (Al) is then replaced by 



3 ira(0i-02)\\ _ 



(A3) 



The magnitude of the fluctuations of M can be characterized by their standard deviation SM 2 = (Af 2 ) — (M) 2 . If 
SM < (M>, then (1/M) ~ 1/ (M) + c5Af 2 / (M) 3 , and Eq. @ becomes 



(A4) 



to be compared with Eq. (Al). If SM <C M and M ^> 1, the correction should be negligible. 

This is why the factor 1/M associated with each particle in the generating function (||) is important. If this factor 
had not been included, the coefficient in front of zz* in C n (z) would be 

(M(M - 1)) ((M* 1 -**)} - (M) ((M) - 1) (e" 1 ^ 1 ) (e~ m ^ 2 ) 

= (M(M - 1)) f/ e ™(^-«\ _ ( e m ^) (e" m ^)) + SM 2 (e" 1 ^ 1 ) (e""*) . (A5) 

In this equation, the second term in the r.-h. s. appears with a coefficient of order SM 2 / (M) 2 with respect to the 
coefficient of the first term. This should be compared with Eq. ( A4 ) , where the coefficient of the second term is of 
magnitude SM 2 / (M) , i.e., much smaller. Thus the weight 1/M in Eq. (J5j) minimizes the effects of a fluctuating 
multiplicity. 

We now turn to the expression of the fourth-order cumulant. For simplicity, we restrict the discussion to a perfect 
detector, and to a fixed multiplicity M. Then, the generating function (G n (z)) is 



(G n (z)) = 1 + ^jp-W 2 (e"^~«) + (M ~ lKil 4 / M3 2)(M " 3) NI 4 (e^+^-W) + 
Inserting this expression in Eq. (^) and expanding to order |z| 4 , one obtains 



s m(0i+</>2— 03 — <t>i) 



M - 1 
4M 



(M — 2)(M — 3) 



,in(4>i+4> 2 - 



, (M-1) 2 
' M 2 



(A6) 



(A7) 



which gives Eq. (Q) in the limit of large M. 

More generally, the cumulants ((exp[m(</>i + ■ ■ ■ + 4>k ~ <f>k+i — ■ • • — (f>k+i)])) derived from the generating equation 
(0) coincide with the quantities of physical interest, that is, only the (k + Z)-particle direct correlations when there is 
no flow, up to subleading terms of relative magnitude 0(1/M). 



APPENDIX B: INTERPOLATION FORMULAS 



In this Appendix, we give interpolation methods to compute numerically the cumulants from their generating 
functions. Generally, one wishes to reconstruct fcj cumulants for integrated flow, i.e., c„{2fc} for k = 1, • • • , fej, and kd 
cumulants for differential flow, i.e., d mn / n {2k + m + 1} for k = 0, • • • , kd — 1. Typical values are ki — 3, m — 1, 2 and 
kd = 2. One thus obtains three independent estimates of integrated flow v n , and two estimates of each v' mn , n . 

The cumulants are defined as coefficients in the power-series expansions of the generating functions C n (z) [Eq. (^)] 
and T> mn / n (z) [Eq. (p7f)]. To extract the cumulants numerically, one first computes the generating functions at the 
points z Pt q — x PA + iy p .q with 

_ ,'2qn 

Xp,q =r y/p COS 
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_ , 2fl7r\ 

y P , q = r \/P sm (Bl) 



for p = 1, • • • , fc max and q = 0, • • • , g max — 1- These equations define a set of fc max x W points in the complex plane, 
where fc max and q max must satisfy fc max > ki, fc^ and g max > 2fcj, 2(fcd + m— 1). With the values of fcj, m, and fc<z given 
above, one can choose for ins tance k max — 3 and q max — 7 or 8. 

The number r in Eqs. ( |Bl| ) is in principle a small number, since we are interested in the behaviour of the generating 
functions near the origin. If ro is too small, however, large numerical errors occur. In practice, the numerical 
simulations presented in Sec. ^ were done with double precision numbers (16 digits), and with the value ro = 1.5. 
A muc h sma ller value was given in |47y . This is due to a rescaling of the variable z by a factor of \fM, as discussed 



in Sec. [V A [see Eq. (|43|)1. In any case, the interpolation should be done with two different values of ro in order to 



check the stability of the results. 

1. Integrated flow 

We denote by C Pi9 the values of the generating function C(z) evaluated at the points (jBlj): 

Cp,q = C n (Zp y q). (B2) 

The cumulants c n {2k} correspond to the terms with k = I in the power-series expansion (g). In order to eliminate 
terms with k ^ I, one averages over the phase of z: 

C p = 2J c p,i- ( B3 ) 

Then, the C p , with p = 1, . . . , ki, are related to the cumulants c„{2fc} with k = 1, . . . , fci by the following linear system 
of fci equations: 



The solution of this system for ki = 3 reads 



C p = E °(^)2 c »i 2fc }> l<P<h. (B4) 



c„{2} = -i ^Ci- ^Ca + iCa 
r o I 2 3 

c«{4} = \ (-5Ci+4C 2 - C 3 ), 

c n {6} = -| (3C 1 -3C 2 + C 3 ). (B5) 



Solving Eqs. (B4) with a larger value of ki provides more accurate values of the first three cumulants, as well as higher 



order cumulants c n {2k}. 

2. Differential flow 

The generating function T> mn j n (z) is complex. From definition (p6[), its real and imaginary parts at the points z Piq 
are 

r M - 9 [D mn/n ^ q) ] = (B6) 

The cumulants d mn / n {2k + m + 1}, defined in Eq. (p9|), are the real parts of the terms proportional to z * k z k + m m 
the power-series expansion (p7[). In order to isolate these terms, one multiplies T> mn / n (z) by z* m , takes the real part, 
and averages over angles: 
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9=0 



9max 



cos ( m ) X p ^q + sin I m I Y PA 



(B7) 



The values of D p for p — 1, . . . , kd are related to the cumulants d mn / n {2k + m + 1} with = 0, . . . , kd — 1 by the 
following linear system of kd equations: 

kd (rn fj5) 2 ( k+m ) 

D P = J2 k {(k + mV d rnn/n{2k + m + l}, 1 < p < k d . (B8) 

k=0 

For kd = 2 and m = 1, the solution of this system is 

rf„/„{4} = \{-2D 1 + D 2 ), (B9) 



while for kd = 2 and m = 2, 



' o 



dan/n{3} = ^(4Di-iD a ) 



d 2 „/«{5} = -T ( ~ 6 Dl + ^ -° 2 ) ■ ( B1 °) 



' o 



APPENDIX C: ACCEPTANCE CORRECTIONS 

We derive here the relations between the cumulants and the flow when the detector has only partial azimuthal 



coverage. As explained in Sec. [II Q , we assume that the class of events used in the flow analysis (usually corresponding 
to a given centrality cut) is selected by means of a detector which has, at least approximately, full azimuthal coverage, 
so that we may assume that the probability distribution of $^ is uniform. 

We wish to recall that the firs t ste p in our accounting for acceptance inhomogeneities is the choice of "nonisotropic" 
cumulants, as explained in Sec. II B: while the simple cumulant deduced from Eq. (||) is valid for a perfect detector, 



a more general definition is that deduced from Eqs. (0) and (||). For "almost" perfect detectors, this first step should 
be enough, and the following results over-sophisticated. 

1. Integrated flow 

Let us describe the detector characteristics (acceptance and efficiency) by a real- valued function A((f>), which is 
the probability that a particle emitted at angle <f> be detected. We choose to normalize this function according to 
J A{4>)d4> = 2tt. The function A(4>) can be expanded into a Fourier series: 

A(4>)= ]T ope**, (CI) 

p— — oo 



where the a p coefficients, Eq. (18), satisfy a_ p = a*, and ao = 1 due to the normalization choice. For a perfect 
detector, a p ^o — 0. 

The distribution of outgoing particles in a given collision is 

+oo 

P(0-$ fl )oc J2 v P e^~*«\ (C2) 

P— — OQ 

where is the reaction plane azimuth of the collision, vq = 1, and the v p = v- p are real valued. 

We now evaluate the generating function of azimuthal correlations (G n (z)}, with G n (z) defined by Eq. (|5|), and 
compute the cumulants as a function of the flow coefficients v p and the acceptance coefficients a p . 
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For a given orientation of the reaction plane, the average value of e for a particle seen in the detector is 



+00 



P( 



d(l>A(<t>)P(4>-$ R ) 



-hoc 



(C3) 



a P v p 



P— — 00 



where n is the harmonic one wants to measure. The denominator is the probability that a pion be detected, which 
depends on if there is flow. 

To obtain the cumulants to leading order in the flow coefficients v p with p 7^ 0, one can linearize this expression in 
v p with p 7^ 0: 



(e- in ^ R ) = a n + J2(a p+n - a n a p )v p e^ R . 

p#0 



(C4) 



In this expression, the term proportional to a n a p comes from the denominator of Eq. (C3). It reflects the <j>^- 
dependence of the probability that a pion be detected. 

The average value (C4) can then be introduced in the generating function (G n (z)). We neglect nonflow correlations 
for simplicity. Then, the angles of the particles are statistically independent and we obtain 



M 



^(a p+ „ - a n a p ) v p e lp * R 

p#0 



M 



e za » +cc exp 



p#0 



dfi dp J Up G 



ip<$> R 



(C5) 



In the second identity, we have assumed that M is large, so that (l + x/M) ~ exp x. Equation (|C5|) can be compared 



to Eq. (|13|), to which it reduces when the acceptance is perfect. When it is not the case, the generating function 
depends in general on harmonics v p with p ^ n. 

The generating function must still be averaged over as in Eq. (|l4|). The first term exp(za„+z*a* ) is independent 
of (fr R and factors out. This term does not contribute to the cumulants of order 2 and higher, since it gives a linear 
contribution to the generating function of the cumulants ([?]), which is In (G n (z)) in the large M limit. 

In the general case, there is no simple analyti c e xpression for the average over 3> R of (G n (z)}. To obtain the 
cumulant at a given order, one must expands Eq. ( |C5| ) to the desired order, and then integrate over <f> R . Then Eq. (j^) 
yields the cumulants. 

Keeping only the first harmonics V\ and v 2 , which correspond to the terms p = ±1,±2 in Eq. (C4), one finally 
obtains for n — 1 (i.e., for a measurement of directed flow) 



Cl {2}= [(1 



2\2 



\a 2 - Oi| 2 ] v\ + [\ai - a 2 a*| 2 + \a 3 - aia 2 \ 2 ] 



(C6a) 



Both harmonics interfere, and one cannot measure v\ and v 2 independently. Similarly, the cumulant to order 4 
becomes 



Cl {4} = - [(1 - M 2 ) 4 + 4 (1 - \ ai \ 2 f\a 2 - a\\ 2 + \a 2 - a 2 | 4 ] v\ 

- [|oi - a 2 a*| 4 + 4|oi — a 2 a*| 2 |a 3 - aia 2 | 2 + |a 3 - aia 2 | 4 ] v 2 . 

For n — 2 (integrated elliptic flow), the first two cumulants are 

c 2 {2} = [|o x - a 2 a^| 2 + |a 3 - aia 2 | 2 ] v\ + [(1 - |a 2 | 2 ) 2 + |a 4 - a 2 2 \ 2 ] v\. 

and 

c 2 {4} = - [\ai - a 2 a{\ 4 -\-A\a\ - a 2 a*| 2 |a 3 - aia 2 | 2 + |a 3 - aia 2 | 4 ] v\ 
[(1 - |a 2 | 2 ) 4 + 4 (1 - \a 2 \ 2 ) 2 \a 4 - a\\ 2 + \a 4 a 2 | 4 ] v A 2 . 



(C6b) 
(C7a) 

(C7b) 



If the acce ptanc e is pe rfect , Eqs. (C6) and (C7) reduce to Eqs. (17a) and (17b), as they should. In the general case, 
Equations (C6a) and (C7a) represent a linear system of equations, which can easily be solved to give v\ and v 2 (or 
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more precisely, vi{2} and «2{2}) as a function of the cumulants. Similarly, Eqs. (C6b) and ( |C7b ) can be solved to 
obtain i>i{4} and ^2{4} as functions of the cumulants. 

We do not include the expressions of the cumulants ci{6} and C2{6}. In addition to the expected terms in vf and 



they involve a term proportional to v\v. 



Note that in going from Eq. (C3) to Eq. (C4), we have kept only the leading order terms in v n , so that the values 
of v\ and i>2 obtained by the above method may have systematic errors of relative order v\ and V2- for instance, if 
the elliptic flow is V2 = 5%, one may find instead 5.25% due to this effect. 

Finally, it should be noted that all corrections for acceptance inhomogeneities, in Eqs. (|C6|) and (|C7|), involve at 



least the squared norm \a p \ of the acceptance Fourier coefficients. Therefore, we believe that if all a. p>0 are smaller 
(in norm) than 0.1, then the procedure exposed in this subsection is superfluous: the relative magnitudes of the 
corresponding acceptance corrections are at most of a few percent, much too small to be significant. Moreover, this 
also applies to differential flow, since as we shall see shortly, in the case of a perfect detector for integrated flow, there 
are no acceptance correction even for v' p . 



2. Differential flow 



For the protons whose differential flow is measured, we introduce two functions A'{ip) and P'(tp ~ ®r) which play 
the same role as A(<fi) and P(cf> — &r) for the pions (please note that A' and P' are NOT the derivatives of A and P). 
Their expansion in Fourier series reads as in Eqs. (CI) and (C2), with a p (resp. v p ) replaced by a' p (resp. v' p ). 

In order to calculate the generating function of the cumulants, Eq. (|26|), we need to evaluate (e jp ^G„(z)), where 
the average is taken over all the detected protons. In computing this average, one must take into account carefully 
that the probability p($_r) that an emitted proton be detected depends on the orientation of the reaction plane &r: 



p{$ R ) = JJ A'^P'^-Qr)^-. 



With this notation, we may write 



(e^G n (z)) 



fa" (e^G n (z)\<Z> R )p(<S> R )d<S> R 



(C8) 



(C9) 



The denominator is equal to 27ra = 2n since we normalize the acceptance function by a' Q = 1. In order to evaluate 
the numerator, we neglect nonflow correlations. Then, the proton and the pions are emitted independently for a fixed 
$>r, so that averages factorize: 



(e^G n (z)\<t> R ) = (e**|$ fl ) (G n («)|*J 



(CIO) 



In this expression, (G„,(z)\§>r) is given by Eq. flC5| ) and (e zp ^|$fl) by an equation similar to (C3), with a p (resp. v p ) 
replaced by a' p (resp. v' p ). Using Eq. (Pq), we may thus write 



„2tt j I +oo 

{^\®r) P {$r) = J ^A'^)P'^)£= J2 <- P <^ R 



(Cll) 



Inserting Eqs. (CIO) and (Cll) in expression (C9), the generating function of the cumulants, Eq. (Eq), reads 



,Cz iq * R {G n {z)\$R)d$R 



V 



p/; 



i (z)= ^ a 'q-p V ' 



^ (G n (z)\<S> R )d<f> 1 



(C12) 



Comparing with Eq. (|31|), to which this equation reduces when the acceptance is perfect, one sees that the cumulants 
involve in general all harmonics v' q . 

Expanding in powers of z and z* and performing the integrals over <&r, one finally obtains for differential directed 
flow 



d V i{2} = » [1 - |ai| 2 + (a' 2 )*(a 2 - af)] v[vi 

+3? K( a i - a 2 a i) + ( a 3)*( fl 3 - 0.10.2)} v' 2 V2, 



(C13a) 
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where to means that one must take the real part: when the acceptance is not perfect, the coefficient is in general 
complex. The higher order cumulant is given by 



di/i{4} = -to 

-to 



"(1 - M 2 ) ((1 - |ai| 2 ) 2 + 2 |o3 - a 2 f) + (a 2 )*(a 2 - a 2 ) (2(1 - K| 2 ) 2 + |a 2 - a 2 | 
a^(aj - a^ai) Ha* - a 2 ai\ 2 + 2 |a 3 - aia 2 | 2 ^) 
+ (a 3 )*(a 3 - aia 2 ) (2 \a\ ~~ a* 2 ai\ 2 + |a 3 - aia 2 | 2 ^ w 2 w| 

Differential elliptic flow measured with respect to integrated directed flow yields 

d 2/1 {3} = K [(1 - K| 2 ) 2 + (oi)*(o2 - a?) 2 ] ^? 

+ 23? [(oi)*(oJ - a*ai)(a 2 - a 2 ) + (a 3 )*(a 3 - a 2 ai)(l - |ai| 2 )] v[viv 2 . 

The differential elliptic flow with respect to integrated elliptic flow is given by 



cZ 2/2 {2} = 3? [l - M 2 + «)*(a 4 - a 2 )] « 2 « 2 

+3? [(a' 1 )*(ai - a 2 a\) + (a 3 )*(a 3 - aia 2 )] v[vi 



and to higher order 
d 2 / 2 {4} 



-to 
-to 



„'„,3 



„'„,3 



(1 - M 2 ) ((1 - M 2 ) 2 + 2 |a 4 - a^| 2 ) + (ai)*(04 - a 2 ) (2(1 - |a 2 | 2 ) 2 + |a 4 - a 2 | 2 ) 
(a' 1 )*(ai - a 2 a^) (|a* - a^ail 2 + 2 |a 3 - aia 2 | 2 ^) 
+ (a 3 )*(a 3 — aia 2 ) (2 |a* — a* 2 a\\ 2 + |a 3 — a\a 2 \ 



„'„,3 



(C13b) 



(C14) 



(C15a) 



(C15b) 



When the acceptance of the detector used for the measurement of integrated flow is perfect, i.e. a p = for p ^ 0, 
these formulas reduce to Eqs. ([mJ) and (|35|), i.e., they do not depend on the differential acceptance coefficients a'. 



APPENDIX D: STATISTICAL ERRORS 



In this appendix, we calculate the statistical fluctuations of the reconstructed integrated and differential flows, 
due to the finite number of events N ovts - More precisely, we calculate the covariance matrices (v n {2k}v n {2l}) — 



(v n {2k}) (v n {2l}) and {v' mn/n {2k + m+ l}v' mn/n {2l + m+l})- {v' mn/n {2k + m + l}) (v' mn/n {2l + m+l}) which 
contain the standard error on each estimate v n {2k} and v' mn / n {2k + m + 1}, and also the linear correlation between 
estimates of different orders. Throughout the appendix, we neglect nonflow correlations and assume that the detector 
is perfect. 

We first introduce some notations. If x is an observable measured in an event (multiplicity, transverse energy, etc.), 
we denote by {x} the average value of x over the available sample of events, which we also call the sampling average: 



{x} = 



1 



Wevts 



(Dl) 



The exact statistical average, corresponding to the limit N cvts — ► 00, will be denoted by (x). Note that ({x}) = (x). 
If y denotes another observable associated with each event, then the covariance of the sampling averages {x} and 

{y} is 



({x}{y})-(x) (y) = 



cvls 



({xy) - (x) (y)), 



(D2) 



where we have used the property that events are statistically independent. When y = x, the square root of the r.-h. s. 
gives the standard deviation of {x} from the exact average (x), which scales like 1 / V iV cv t s . This result will prove 
useful later on. 
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1. Two-point correlation function 



As explained in Sec. [IB, one constructs for each event a generating function G n (z), defined by Eq. (|5|). The 
generating function of azimuthal correlations is the statistical average of this function, (G„(z)}. This quantity, given 
by Eq. ((L5|), does not depend on the phase of z because of isotropy. On the other hand, one measures experimentally a 
sampling average {G n (z)}, which generally has a (weak) dependence on the phase of z, due to statistical fluctuations. 

In the following, we shall need to evaluate the statistical fluctuations of the sampling average {G n (z)} around the 
true statistical average (G n (z)). These statistical fluctuations are characterized by the two-point correlation function 



({G n (z)} {G n (z')}) (G n (z)) (G n (z>)) = - J- ((G n (z)G n (z>)) - (G n (z)) (G n (z>))) , 

JVovts 



(D3) 



where we have used Eq. (D2). To evaluate the r.-h. s. of Eq. (D3), we first perform the average for a fixed orientation 
of the reaction plane Using Eqs. (^) and (|l2|), we obtain 



(G n (z)G n (z')\$ R ) = 1 + 



(z + z')v n e~ m * R + c.c. zV + c.c. 



M M 2 
~ exp \^z + z')v n e - in *« + ^+ c.c. 

where c.c. denotes the complex conjugate. One must then average over $f>: 



M 



{G n (z)G n (z')) = 



2tt 



(G n (z)G n (z')\$ R ) 



(D4) 



(D5) 



This function is invariant under a global rotation (z,z') — > (ze m , z'e ia ), but depends on the relative phase e l6 
z*z' j\zz'\ between z' and z. We can expand it in Fourier series with respect to 9, in the form 



(G n (z)G n (z')) = S m (\z\,\A) 



I \ m, 



\zz'\ 



(D6) 



where the Fourier coefficients are given by 



2jv 



(10 



g m {\z\ t \A)= I ^e- im<> (G n (\z\e ie )G n (\z'\)) 



d9 

o 2^' 



d<f>! 



■ exp 



(\z\v n e *(«-»*«) + \z'\v n e- in ** + J^l, 



(D7) 



If flow is small, more precisely if v n <C 1/s/M, we can set v n = in this equation. Then, the integral over $^ is 
trivial, while the integral over 9 yields 



Gm(\z\, \z'\) = I„ 



2\zz'\ 
M 



When flow is larger, the integrations can be performed using the following identity, valid for real x: 



exp(xe^+xe"^) = ^ I q (2x) . 



(D8) 



(D9) 



q— — oo 



Note that 7_ g (2x) = I q {2x). Applying Eq. (DE) to the three terms in the exponential in Eq. (D7), one obtains 

'2\zz'\ 



Q m {\z\,\z'\) = I q (2\z\Vn)I q (2\z'\v n )I t 



q— — oo 



M 



(D10) 



When v n — 0, all terms in the sum vanish but q = 0, and one recovers Eq. (D8). 
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2. Integrated flow 



Experimentally, the generating function of the cumulants C n (z) is obtained by replacing (G„(z)) with the sampling 
average {G n (z)} in Eq. (R). The cumulants c n {2k} used to estimate the integrated flow are then obtained through a 
power-series expansion. We first introduce the notation 



Qk Ql 



2 = 



With this notation, the cumulant c n {2k}, defined by Eqs. (|§|) and can be written 

c„{2fc} = C n (z)\ zkz , k . 

The statistical fluctuations of the cumulants c„{2fc} is characterized by their covariance matrix: 
(c n {2k}c n {2l}) - (c n {2k}) (c n {2l}) = (C n (z)C n (z')) - (C n (z)) {C n {z'))\^ 2k]zl?l 



(Dll) 



(D12) 



(D13) 



Let us evaluate the r.-h. s of this equation. We assume that the multiplicity is large enough, so that C n (z), defined 
through Eq. (Q), reduces to ln{G„(z)} for the sample of events. Expanding the logarithm to first order around 
(G n (z)), one obtains: 



(C n {z)C n (z'))-(C n {z)) (C n (z')) = 



({G n (z)} {G n (z>)}) - (G n (z)) (G n (z>)) 

(G n (z)) (G n (z>)) 
(G n {z)G n (z')) 



N cvts \{G n (z)) (G n (z')) 



(D14) 



where we have used Eq. (D3). According to Eq. ( D13 ), we must isolate the terms proportional to |z| 2fc |z'| 2i in the 
power series expansion, i.e., terms which do not depend on the phases of z and z' . The denominator, given by 
Eq. (|l5|), depends only on \z\ and \z'\. On the other hand, the numerator depends on the relative phase of z and z' . 
The isotropic part is the term m = in the Fourier expansion (D6). We thus obtain: 



(c„{2fc}c„{2Z}) - (c„{2fc}) (c„{2/}) 



Go(\z\,\z'\) 



N evts \I (2\z\v n )I {2\z'\v n ) 



(D15) 



|z| 2fc |z'| 2 



where <?o(M, \z'\) is given by Eq. ( D10| ) with m = 0. 

The estimate of integrated flow v n {2k} is obtained by expanding Eq. jl^ ) to order \z\ 2k . Using Eq. ( D12 ), it is 
related to c n {2k} by 



c n {2k} = m/ (2|*K{2fc})L „ ,„ 



(D16) 



Using this equation, one easily relates the covariance matrix of the estimates v n {2k} to that of the corresponding 
cumulants. Noting that 



ln/ (2|2|K + Sv)) - In J (2|zK) = 2 \z\ 



Io(2\z\v r , 



(D17) 



we obtain from Eq. ( D16 ) 



(c n {2k}c n {2l}) {c n {2k}) {c n {2l}) = 4\zz'\ Il{2lzlVn) h{2 \ z '\ v ^ 



'/o(2|«|i;„) Io(2|*'K) 



({v n {2k}v n {2l}) - (v n {2k}) {v n {2l})) . 



|z| 2fe |z'| 2 



(D18) 



Using the expression of the covariance matrix of the cumulants obtained above, Eq. ( D15| ), one thus obtains a compact 
expression for the covariance matrix of the estimates v n {2k}. 

Before giving explicit results for the lo wes t order cumulants, let us dis cuss the weak flow and strong flow limits. If 
v n -C 1/%/M, <?o(Mj \ z '\) reduces to Eq. (|D8|), and the r.-h. s of Eq. ( DIE ) depends only on \zz'\, thus terms with k ^ I 
vanish: correlations between different cum ulants vanish in this limit, so that different estimates of v n are uncorrelated. 
If v„ » 1/VM, one can expand Eq. (DIP) to order 1/M, which yields 
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2 \zz'\ 

g (\z\, \z'\) = J (2|zK)/ (2|z'K) + -^h(2\z\v n )h{2\z'\v n ). (D19) 



Inserting this expression into Eq. (D15), and comparing with Eq. ( Dig ), we obtain 

(v n {2k}v n {2l}) - (v n {2k}) (v n {2l}) = 1 . (D20) 

zivi Jv ovts 

In this limit, all estimates v n {2k} coincide and the error on the integrated flow is l/^/2MN cvts , independent of A:. This 
result can be easily understood: when v n is large compared to the reaction plane 3>r can be reconstructed with 

very good accuracy. Then, the integrated flow v n {2k} can be obtained as the average over all particles of all events of 
cos[n(4>j — $r)]- Since the total number of particles is MN evts , one evaluates (l/MN cvts ) X) J= i cvts cos[n(<f)j — $#)], 
and the average value of the square of this quantity is 1/(2M N cvts ) for random angles. Thus the resulting statistical 
error is 1/ y / 2MN cvt s- 

Finally, in more general case when v n and 1/yM are of the same order of magnitude, explicit expressions for the 
lowest order estimates are obtained by expanding Eqs. (D15) and ( D18| ) in power series of \z\ 2 and \z'\ 2 . The standard 
deviations on the first order estimates are given by 

(Sv n {2}) 2 ee {v n {2} 2 ) - (v n {2}) 2 = ' ! " 2 V 



2MN cvts 2 X 

(6v n {4}f ee < v „{4} 2 ) - (v n {4}f = 



2_ / r,!2\ /.. r^xa _ 1 l + 4 X 2 +x 4 + 2 X 6 



2MN cvts 2 X e 

fx saw 2 -/ /«ia\ / srw 2 1 3 + 18 X 2 + 9x 4 + 28 X 6 + 12 X 8 + 24 x 10 . 
(Sv n {6}) = (v n {6} ) - W8}) - 2^35 (D21) 



with x 2 = Mi^. In the limit x> 1, the y resu lts reduce to Eq. ( D20 ). 

The values of Sv n {2k} given by Eqs. (D21) are plotted in Fig. |2j as a function of x- We have taken iV cvts = 10 s 
events, with multiplicity M = 200 each. With these values, \ — 1 corresponds to v n ~ 7%. Some comments on these 
results: 

• The statistical error increases rapidly as v n decreases: if v n = 3% (i.e., x = 0.42), the statistical error using 
the fourth-order cumulant is reasonably small, less than 0.2%. But if v n — 1.5% (x = 0.21) the statistical error 
becomes of the same magnitude as the flow itself. Such small values of the flow can be studied only with a 
higher multiplicity and/or a very large number of events. 

• While the lowest order estimate v n {2} always has the smallest statistical error, the statistical error on u„{6} is 
slightly smaller (by at most 5%) than the error on i> n {4} for \ > 1- 

The correlation between estimates of different orders is given by 
(«„{2K{4}) - (v n {2}) (v n {4}) = — 1 

Z1V1 iV evts 

(v n {2}v n {6}) - (v n {2}) (v n {6}) = — 1 

INI iVevts 

(v n {4}v n {6}) - K{4}) (v n {6}) = — 1— + 9 * 6 8 +2x8 (D22) 

2MN evts 2% s 

Figure |2| displays the linear correlation C{ 2 fc},{2;} between v n {2k} and v n {2l} with k ^ I. We recall that the linear 
correlation c between two random variables x and y is defined as 

(*y> -<*><») (D23) 



y/(x*) - (x) V(y 2 > - (y) 2 

c always lies between —1 and 1; these two limiting cases corresponding to a linear relation y = ax, while c = if x 
and y are uncorrelated. As expected from the discussion following Eq. (D20), different estimates are uncorrelated if 
X <C 1, but the correlation becomes stronger and stronger as the resolution x increases. 



2G 



0.5 



1.5 



2.5 



FIG. 2. Top: statistical errors on v n {2}, v n {A} and v n {6} for iV cv ts = 10 5 and M = 200, as a function of v = 
Bottom: linear correlation C{2*},{2!} between each pair of estimates v„{2k} and v„{2l}, with 1 < k < 3 and 1 < Z < 



Finally, we would like to point out that the statistical fluctuations of the estimate v n {2k} around the true value 
v n are not Gaussian. Indeed, it can be shown that the fluctuations of the cumulants are generally Gaussian. Thus, 
according to Eqs. (|l7|), the fluctuations of the variable £ = v n {2k} 2k are Gaussian, but not the fluctuations of v n {2k} 
itself. We may write 



dN 



1 



: exp 



(d-vl k ) 2 
2a 2 



where a is the standard deviation on £. It is related to the deviation Sv n {2k}, Eqs. (D21), by 

l Sv n {2k}. 



2kvi 



(D24) 



(D25) 



In order to illustrate the non-Gaussian character of the fluctua tions of v n {2k}, we display in Fig. ^ the distribution 
of v n {2}, w„{4} and v n {6} when the a verage value of £ in Eq. (D25) is only one standard deviation above zero, i.e., 



when a = w„ . Then, using Eq. ( D24 ), the probability that £ < is about 16%, in which case v n {2k} is undefined. 



One notes that for the 84% remaining cases, the distribution of v n {2k} becomes narrower as 2k increases. 



2.5 



1.5 



0.5 





"n{2) 








/ i ''n{6} - - - - 







0.5 



1.5 



"n(2k}/v n 

FIG. 3. Distribution of v n {2}, u n {4} and i; n {6}, scaled by the theoretical value v n - 
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3. Differential flow 



The generating function of the cumulants used for differential flow, D p / n (z), is the ratio of two quantities: the 
numerator of the r.-h. s. of Eq. ( |26| ) is evaluated from a sample of AT' protons, while the denominator is calculated 
from the same sample of events as for the integrated flow, that is, N cvts . Since the measurement of differential flow is 
usually performed in a narrow phase-space window, we neglect the contribution of the denominator to the statistical 
error. We thus write 



{e^G n (z)} 

(Gn(z)) ' 



(D26) 



where the denominator is given by Eq. ( |l5| ) 
The cumulants d. 



mn/n{2k + m + 1} defined by Eqs. ( p7[ ) and ( |29| ) can be written in the form 

1 



Cn/«{2fc + m + 1} = - (V mn/n (z) + V* mn/n { 



z*) 



(D27) 



z k + m z *k 



We now evaluate their covariance matrix. For simplicity, we shall assume that there is no differential flow, i.e., that 
t he ij) distribution is isotropic, so that CD p / n (z)) — 0, and only fluctuations remain. Using Eq. (D2) and the definition 
( |D2q ), we thus obtain 



{V p/n {z)V p/n {z')) = 

(v p/n {z)v; /n {z')) 



1 (G n (z)G n (z')) 



N' (G n (z))(G n (z>))' 



(D28) 



We recogni ze in the last equation the two-point function studied in Appendix D 1 . The covariance matrix of the 
cumulants ( D27 ) can thus be written 



/, /9 . , m , n . r 97 . m ■ n \ 1 ( i G n(z)Gn(z'*)) , (G n (z*)G n (z')) 

{dmn/nV-K + TO + 1 }d mn/n {2l + TO + 1}) = — 



47V' V (G n {z)) (G n (z>*)) (G n {z*)) (G n (z>)) 
1 (G n (z*)G n {z')) 



z k + m z *k z ll + -m. z i*l 



z k + mz *k z ll + m z l*l 



2N> (G n (z)) (G„(z')> 

Using Eq. ( p^| ) and the Fourier decomposition of the two-point function (D6), we obtain 

(d mn/n {2k + m + \}d mn/n {2l + to + 1}) = -L J^l^ ,, . { , 

2AT' Jo(2|^|Vn)io(2|r|Un) \\zz '\, 



f zz' 


^ in 






[\zz'\ 




z k + 


m z *k z H + m z f* 



(D29) 



(D30) 



where ^ TO (|z|, \z'\) is given by Eq. (pHj). 

The estimate of differential flow v' m ^ n {2k + m + 1} is obtained by expanding Eq. ( |3l| ) to order z k + m z * k _ Using 
Eq. (D27), it is related to d mn / n {2k + to + 1} by 



t {2k + to + f } = y^pr^y f 7~T 
I Q (2\z\v n ) \\z\ 



,k + m z *k 



Using this equation, one easily relates the covariance matrix of the estimates v' mn ^ n {2k 
corresponding cumulants, given by Eq. (D30): 

I m {2\z\v n ) I m (2\z'\v n ) 



(D31) 

1} to that of the 



J (2|zK) I (2\z'\v n ) 



\z\ 2k \z'\ 



(<„/„{2fc + m + l}v' mn/n {2l + m+ 1}) = ^ - 



Gm(\z\,\z'\) 



{2\z\v n )I (2\z'\v n ) 



\z\ 2k \z'\ 2 ' 

(D32) 



Let us discuss the weak flow and strong flow limits. When v n <C 1/yM, Q m (\z\, \z'\) is given by Eq. (|D8|), and 
terms with k ^ I vanish: as in the case of integ rated flow, correlations between estimates of different orders vanish in 
this limit. When v n ^> 1/VM, expanding Eq. ( D10| ) to leading order in 1/M, we obtain 



Qra ( | Z 

I, \ z '\) = Im(2\z\v n )I m (2\z / \v n ). 



(D33) 



2cS 



Then, the covariance matrix, Eq. (D32), reduces to 

(<m/n{2* + m + l}v' mn/n {2l + m + 1}) - (v' mn/n {2k + m+ 1}) (v' mn/n {2l + m+ 1}) = JL. (D34) 

This result can be simply understood by repeating the argument used in the case of integrated flow. 

In the more general case whe n v n and 1 / ^/M are of the same order of magnitude, the following lowest order formulas 
are derived by expanding Eq. (D32) in power series. For m = 1: 

1 I + x 2 



;</„{2} 2 ) - (</n{2>) 

(</„{2K /n {4}) - (< /n {2}) (< /n {4} 

;< /n {4} 2 ) - (</„{4}) 



27V' x 2 
1 

27V 7 

1 2 + 6 X 2 + X 4 + X 6 
2AT' x 6 



(D35) 



with x 2 = Mvff For m = 2: 



; 2n/; 



„{3}^ 



; 2n/n{ 3 }/ — 2^y/ 



1 2 + 4 X 2 + X 4 



X 

U 2n/n{ 3 Kn/n{5}) ~ ( U 2n/«{ 3 }) ( V 2n/n{ 5 }) = ^7 ~~ 

2 



«2„/ J5} 2 ) - (»J n/n {5})' 



1 6 + 24x 2 + 9x 4 + 10x 6 + 4x 8 
27V 7 



4 X s 



(D36) 



When x ^S> 1, these results reduce to (D34). Figure |4| displays the variation with \ of the standard errors on the 
various estimates and the linear correlation between cumulants of different orders. The behaviour is qualitatively the 
same as for integrated flow. 



0.1 





\ \ 8v 2n/n{ 3 } 




\ \ Sv 2n/n{ 5 } 






c {2},{4} 





0.5 1 1.5 2 0.5 1 1.5 2 2.5 

X 

FIG. 4. Left (resp. right): statistical properties of the estimates of v' n / n (resp. v' 2n , n ). Top: statistical errors on v' n / n {2} and 

u n/n{4} (resp. «2n/n{3} and v' 2 „/ n {5}) f° r N' = 5 x 10 5 and M = 200, as a function of \ = VnVM. Bottom: linear correlation 
between v' n/n {2] and < /n {4} (resp. « 2n/n {3} and « 2n/ „{5}). 
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